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Summary 

An analysis method based on a deformation (as opposed to damage) approach has been developed to 
model the strain rate dependent, nonlinear deformation of woven ceramic matrix composites with a plain 
weave fiber architecture. In the developed model, the differences in the tension and compression response 
have also been considered. State variable based viscoplastic equations originally developed for metals 
have been modified to analyze the ceramic matrix composites. To account for the tension/compression 
asymmetry in the material, the effective stress and effective inelastic strain definitions have been 
modified. The equations have also been modified to account for the fact that in an orthotropic composite 
the in-plane shear stiffness is independent of the stiffness in the normal directions. The developed 
equations have been implemented into a commercially available transient dynamic finite element code, 
LS-DYNA, through the use of user defined subroutines (UMATs). The tensile, compressive, and shear 
deformation of a representative plain weave woven ceramic matrix composite are computed and 
compared to experimental results. The computed values correlate well to the experimental data, 
demonstrating the ability of the model to accurately compute the deformation response of woven ceramic 
matrix composites. 


Introduction 

Ceramic matrix composites are gaining wider usage in the aerospace industry due to their ability to 
withstand high temperatures while still providing relatively high stiffness and strength levels as compared 
to monolithic ceramics. Woven ceramic matrix composites with a 2-D plane weave fiber architecture are 
commonly used due to their ability to provide high levels of stiffness and strength in multiple directions. 
One area of interest to NASA is studying how these materials respond to high strain rate impact loads. 
However, these ceramic composites have several characteristics in the macroscopic deformation response 
which make modeling the deformation response under impact conditions challenging. The stress-strain 
curve for these materials are nonlinear, particularly in tension and shear, due to factors such as matrix 
microcracking (refs. 1 and 2), and the deformation response changes with strain rate (refs. 3 and 4). 
Furthermore, the uniaxial tension and compression responses are different (refs. 1 and 5), which implies 
that on the macroscopic scale the hydrostatic stresses have an effect on the material response. 

Previous efforts have focused on using damage mechanics approaches to model the nonlinear 
deformation of woven ceramic matrix composites (refs. 1,5, and 6). However, there are several features of 
these models which indicate that examining alternative methods of simulating the material response 
would be beneficial. First of all, in these previous models strain rate dependence of the material behavior 
is not explicitly accounted for within the model. However, previous studies on these materials 
(refs. 3 and 4) have indicated that these materials indeed have a nontrivial rate dependence. Furthermore, 
in the previous studies referenced above all of the nonlinearity of the deformation response is assumed to 
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be due to damage mechanisms instead of deformation mechanisms. Providing an alternative modeling 
capability in which the nonlinearity in the deformation response is assumed to be due to deformation 
mechanisms could allow for a decoupling of the deformation and damage/failure response, which could 
facilitate improved modeling of both the deformation response and the prediction of ultimate failure. 

At NASA Glenn Research Center, efforts have been underway for several years to model the strain 
rate dependent, nonlinear deformation of polymers with the effects of hydrostatic stresses included. In 
these efforts, the Bodner-Partom viscoplasticity model based on the state variable approach (ref. 7), 
originally developed for metals, has been modified in order to incorporate the hydrostatic stress effects 
known to be present in polymers (ref. 8). These variations include modifying the effective stress and 
effective inelastic strain definitions, based on an application of the Drucker-Prager yield criterion (ref. 9), 
to account for the effects of the hydrostatic stresses. 

In this paper, the constitutive equations with an associative flow rule mentioned above have been 
applied and modified in order to model the nonlinear, strain rate dependent, deformation of plain weave 
woven ceramic matrix composites, in which the differences in the tensile and compressive response are 
accounted for within the equations. Due to the plain weave fabric geometry of the materials, for the 
purposes of this study the material is assumed to be isotropic in the normal directions, even though these 
materials are composites and not homogenous materials. However, since the woven composites are 
indeed orthotropic materials on the macroscopic level, the shear response is assumed to be uncoupled 
from and independent of the normal response. The full theoretical development of the constitutive 
equations will be described in this paper, along with the procedures to be used to determine the material 
constants. The procedures used to implement the equations within a commercially available transient 
dynamic finite element code, LS-DYNA (ref. 10), through the use of user defined subroutines will also be 
described. Finally, the tensile, compressive, and shear deformation of a representative plain weave woven 
ceramic matrix composite will be computed and compared to experimental results to demonstrate the 
ability of the model to capture the features of the deformation response of plain weave woven ceramic 
matrix composite materials. 


Theoretical Development 

To analyze the nonlinear, strain rate dependent deformation of plain weave woven ceramic matrix 
composites, the Bodner-Partom viscoplastic state variable model (ref. 7), which was originally developed 
to analyze the viscoplastic deformation of metals above one -half of the melting temperature and has been 
successfully used in such instances, has been modified. This particular model was chosen due to previous 
success in using a variation of this model to analyze the strain rate dependent, nonlinear deformation of 
polymer materials (ref. 8). In state variable models, a single unified strain variable is defined to represent 
all inelastic strains (ref. 11). Furthermore, in the state variable approach there is no defined yield stress. 
Inelastic strains are assumed to be present at all values of stress, the inelastic strains are just assumed to 
be very small in the “elastic” range of deformation. State variables, which evolve with stress and inelastic 
strain, are defined to represent the average effects of the deformation mechanisms. 

An important point to note is that frequently, in state variable models, each of the state variables and 
material constants have a fairly explicit link to specific deformation mechanisms, such as dislocation 
pileups in metals (ref. 11). Here, due to the complex nature of the material deformation mechanisms in 
woven ceramic matrix composites the linkage between the state variables and material constants to 
specific deformation mechanisms is not nearly as well defined. For this work, an internal stress state 
variable is defined to represent phenomelogically the resistance to inelastic deformation, a state variable 
is defined to represent the effects of hydrostatic stresses and a state variable is defined in order to scale 
the shear stresses to ensure a consistent effective stress. The material constants are characterized in order 
to correlate with the experimentally observed deformation of the material, and are not tightly linked to 
specific deformation mechanisms. However, the equations developed here employ a fairly simple 
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formulation with a minimum of material constants that are fairly simple to characterize, and capture the 
main features of the defoimation response. 

For this study, temperature and moisture effects have been neglected, as only room temperature data 
are currently available. All of the nonlinearity and strain rate dependence observed in the material is 
assumed in this work to be due to inelastic defoimation, where in reality the nonlinearity is most likely 
due to a mixture of defoimation and damage. Ultimate damage and failure criteria will need to be added 
to the model at a later date. Small strain theory is assumed to apply, and phenomena such as creep, 
relaxation and high cycle fatigue are not accounted for within the equations. Furthermore, as mentioned 
previously, in reality plain weave woven ceramic matrix composites are orthotropic materials, as is 
typical of composites. Flowever, due to the plain weave architecture of the fibers, the material can be 
assumed to be isotropic in the normal directions, with independent shear responses. 


Associated Flow Equation and Evolution Equations 

In order to derive the associated flow equation for the components of the inelastic strain rate tensor 
for plain weave woven ceramic matrix composites with hydrostatic stress effects (different responses in 
tension and compression) included, a procedure similar to that employed to derive the Prandtl-Reuss 
equations in classical plasticity is utilized (refs. 9 and 11). The formulation employed in this study is 
based on that used by Pan and co-workers (refs. 12 and 13) and FIsu, et al (ref. 14). First, an inelastic 
potential function based on the Drucker-Prager yield criterion (ref. 9) is defined as 


/ = 



+ aa kk 


( 1 ) 


where o kk is the sum of the normal stress components (assuming the usual rules of indicial notation 
(ref. 11)) and is equal to three times the hydrostatic stress. The variable a is a state variable which 
controls the level of the hydrostatic stress effects. The term “<xc kk ” incoiporates the effects of hydrostatic 
stresses into the inelastic potential function. The variable J 2 is a modified version of the second invariant 
of the deviatoric stress tensor. The modification is required due to the fact that for the orthotropic 
composites examined here, the shear terms are independent of the normal terms in the stiffness matrix. 
Therefore, in order to ensure a consistent effective stress (meaning a constant effective stress vs. effective 
strain curve independent of the load type), the shear stresses must be multiplied by a scaling state variable 
p. The modified version of the deviatoric second invariant is given as follows 


* 1 

J 2 = — (<T i -C 22 Y +(^22 -c 33)’ + ( a 33 -(7 ll)“ +P a f2 + °13 +a 23 
6 


( 2 ) 


where Gy are the components of the stress tensor. Note that for the current work all three shear stress 
components are scaled by the same value of p. In actuality, each shear stress component may require a 
unique scaling factor. However, for the current work only in-plane loads are considered, and in the finite 
element implementation the model is applied using thin shell elements only at the current time, which 
results in the transverse shear stresses being approximate at best. Furthermore, at the current time only in- 
plane shear data is available to be examined. As a result, only the correct scaling of the in-plane shear 
stresses is considered significant for this work. The effects of possible errors in the scaling of the 
transverse shear stresses are not considered at this time. The full effects of the transverse shear stresses 
will be considered in future work. 

The components of the inelastic strain rate tensor, £ , are assumed to be proportional to the 

u 

derivative of the potential function with respect to the components of the stress tensor, Gy, as follows 
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where A, is a scalar rate variable. By taking the specified derivative of the potential function, the result 
becomes 


y s l 

dOij 


+ a.5,; 

y 


( 4 ) 


where S, ; is the Kronecker delta and Sj are the modified components of the deviatoric stress tensor. For 
the normal stress components, S tj is equal to Sy, the components of the deviatoric stress tensor. Flowever, 
to maintain consistency, for the shear stress components the terms are scaled as follows 

S*j=(Shj ( 5 ) 

where (3 is the shear scaling state variable discussed earlier. 

The next step in defining the flow rule is to define the effective stress and the effective inelastic strain 
rate. The effective stress, o e , is defined as follows 

o e =Sf = + So.a kk . (6) 


To determine the value of the scalar A in terms of the effective inelastic strain rate, , the principal of 
the equivalence of the inelastic work rate is employed. With the aid of equations (3) and (4), the inelastic 
work rate, W 1 , can be expressed by the following 


W 1 =v A =°ij 


S* 


, — + a8,7 

2-JA ] 


A = a e i e . 


(7) 


After applying the effective stress definition given in equation (6) and simplifying, it can be shown that 


A = Vse; 


(8) 


By substituting equation (8) and equation (4) into equation (3), a definition for the components of the 
inelastic strain rate tensor can be given as 


■ I 


i = A>i\ 


f * 

S» 

J +a8 


v 


2-Ah 


v 


(9) 


By utilizing equation (7), it can be shown that the effective inelastic strain rate can be defined as 
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where el is the effective deviatoric inelastic strain rate and £ l m is the mean inelastic strain rate. The 
inelastic shear strain rates need to be scaled by a constant scale factor k due to the fact that for an 
orthotropic material the shear terms in the material response are independent of the normal terms. When k 
equals one, which would happen for isotropic materials, the inelastic strain rate definition given here 
matches the effective inelastic strain rate definition given by Pan and co-workers (refs. 12 and 13). 

By following a similar format to the Bodner-Partom model (refs. 7 and 11), the following definition is 
specified 


V3 ./ _ 1 

— e e =D 0 exp -- 


' z ^ 2 " 


V c ey 


( 11 ) 


where D 0 and n are material constants. D„ represents the maximum inelastic strain rate, and n controls the 
rate dependence of the material. Z is an isotropic state variable which represents the resistance to inelastic 
deformation (internal stress). By substituting equation (11) into equation (9), the final form of the flow 
equation is determined to be 


4 = 2 D 0 exp| 


1 

' Z ' 

2 n 

f s a 1 

i a /S 

2 



■ i 1 XU 77 

(NT ') 


( 12 ) 


Note that the elastic components of the strain rate must be added to the inelastic strain rate to obtain the 
total strain rate. 

The rate of evolution of the internal stress state variable Z, the hydrostatic stress effect state variable 
a and the shear scaling factor state variable (3 are defined by the equations 


Z = q(Z\ -Z)e e 

(13) 

a = g(oq -a)el 

(14) 

P = ^(Pi-P)4 

(15) 


where q is a material constant representing the “hardening” rate, and Z b oq, and (3 1 are material constants 
representing the maximum values of Z, a, and (3 respectively. The initial values of Z, oq and (3 are defined 

by the material constants Z„, oq and (3 0 . The term e l e in equations (13), (14), and (15) represents the 
effective deviatoric inelastic strain rate, which was defined in equation (10). An important point to note is 
that in the original Bodner model (ref. 7), the inelastic work rate was used instead of the effective inelastic 
strain rate in the evolution equation for the internal stress state variable. However, for this work the 
inelastic strain rate was deemed easier to work with from both computational and characterization points 
of view, particularly in the incorporation of hydrostatic stress effects. Furthermore, the effective inelastic 
strain rate has been used in other state variable constitutive models (ref. 1 1). Since hydrostatic stress 


NASA/TM— 2004-213125 


5 



effects were not considered in the original Bodner model (ref. 7), the evolution equation for a is new to 
this work. Furthermore, since these equations were originally developed for homogenous isotropic 
materials and not composites, the evolution equation for (3 is also new to this work. The state variables 
a and (3 are assumed to evolve in the same manner as the state variable Z. By using this assumption the 
value of q used in equations (14) and (15) can be shown to be the same as the value of q used in 
equation (13). 


Determination of Material Constants 

The material constants that need to be determined include D 0 , n, Z„, Z h oq, oq, (3 0 , (3 1 , k, and q. These 
values, along with the elastic modulus, shear moduli and Poisson’s ratio are the input parameters required 
for the model. The procedure to be used is summarized here. More details on the procedure as it was used 
for isotropic polymers can be found in Goldberg, et al. (ref. 8). The constants are determined using tensile 
and compressive stress-strain curves obtained from a set of constant strain rate tests, along with a shear 
stress-shear strain curve obtained for at least one strain rate. Each tension and compression test is 
conducted at a different total strain rate. 

The first step in the process is to obtain the values of oq and a.„. To accomplish this task, equation (6) 
is used in combination with stress-strain data from constant strain rate uniaxial tension and compression 
tests. The values of these constants are assumed to be rate independent, so the results from only one strain 
rate need to be used to find the needed parameters. In practical application of the methodology, the 
uniaxial tension and compression tests used do not have to be at the exact same effective strain rate. As 
long as the effective strain rates from the two tests are approximately equal, the values obtained have been 
found to be valid, particularly if quasi-static strain rates are examined. Two stress values from the tensile 
and compressive stress-strain curves are required in order to find these two constants. First of all, the 
“saturation stress”, the stress level at which the stress-strain curve flattens out and becomes horizontal, is 
required. At this stress level, the inelastic strain rate is assumed to be equal to the total applied strain rate. 
If the actual stress-strain curve terminates before the saturation stress is reached (as often happens in 
compression), the curve must be extrapolated to obtain an approximate value of the saturation stress. 
Secondly, the stress level at which the tension and compression curves becomes nonlinear (the “yield” 
stress), is required. One important point needs to be considered when determining the stress level 
designating the onset of nonlinearity. Frequently, in woven ceramic matrix composites, the tensile 
stiffness across nearly the entire stress-strain curve tends to be significantly softer than the compressive 
stiffness, most likely due to the fact that in tension significant matrix microcracks develop almost 
immediately upon loading, which weakens the material. As a result, in cases where this phenomena 
occurs the initial portion of the compression curve is assumed to represent the “elastic” regime of the 
composite material, and in tension the onset of nonlinearity in the stress-strain curve is assumed to occur 
at a very small, almost negligible stress value. 

Given these definitions, equation (6) can now be used to determine the values of oq and oq The 
primary assumption used at this point (and assumed implicitly in equation (6)) is that the effective stress 
at saturation under uniaxial tensile loading at a particular strain rate is equal to the effective stress at 
saturation under compression loading at the same equivalent strain rate. Likewise, the effective stress at 
the point the stress-strain curve becomes nonlinear under tensile loading is equal to the effective stress at 
the point the stress-strain curve becomes nonlinear under compression loading. Therefore, assuming the 
value of a at saturation is equal to oq, and the value of a at the point the stress-strain curve becomes 
nonlinear is equal to oq, the following equations are obtained for the case of having data from uniaxial 
tension and compression tests 


’ st 


(l + V3oq )= a sc (l - V3oq ) 


(16) 
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where G„ and G sc are the tensile and compressive stresses at saturation, respectively, and G„/, and G„ fc are 
the tensile and compressive stresses at the point where the respective stress-strain curves become 
nonlinear. Therefore, by substituting the tensile and compressive stresses at saturation into this equation, 
followed by the tensile and compressive stresses at the point where the respective stress-strain curves 
become nonlinear, the required constants can be determined. At this point, by applying equation (6) for 
the tensile and shear curves, using the appropriate “saturation” and “yield” stresses in tension and shear, 
the values of (3 0 and (3! can be determined from the following equations 


O s ,(l + V3(Xi)= 

(18) 

^ nil (l F V-3C( 0 | — yj 3(3 0 X n i 

(19) 


where o st and x s are the tensile and shear stresses at saturation, respectively, and G„* and x„/ are the tensile 
and shear stresses at the point where the respective stress-strain curves become nonlinear. 

The values of D 0 , n, and Z x are characterized as follows. The value of D 0 is currently assumed to be 
equal to a value of 10 4 times the maximum applied strain rate, which correlates with the maximum 
inelastic strain rate. To determine the values of n and Z,, first equation (12) is simplified to the case of 
uniaxial tensile loading, leading to the following expression 


t 1 = 2 Dq exp 


\ 2/7 


(1 + V3a)G 


V3|< 


- + a 


( 20 ) 


where z 1 is the uniaxial tensile inelastic strain rate in a constant strain rate tensile test, o is the uniaxial 
tensile stress, and the remainder of the terms are as defined earlier. The case of uniaxial tensile loading is 
used to determine these constants since for woven ceramic matrix composites the tensile stress-strain 
curves obtained experimentally are more likely to display a defined “saturation” stress than the 
compression curves, which as shown below is crucial for determining the material constants. 
Furthermore, equation (20) can be rearranged as follows 


- 2 In 


2 Dr 


s 


+ a 


, 2 n 


(1 + V3oc)g 


2) 


( 21 ) 


The required constants are determined from a set of tensile stress-strain curves obtained from constant 
strain rate tests. Each curve in this set is obtained at a different constant strain rate. Data pairs of the total 
strain rate (equal to the inelastic strain rate at saturation) and saturation stress values from each curve are 
taken. For each strain rate, the data values are substituted into equation (21), after the natural logarithm of 
both sides of the equation is taken, and represent a point on a master curve. The number of points in the 
master curve equals the number of strain rates at which the tensile tests were conducted. A least squares 
regression analysis is then performed on the master curve. As can be shown from equation (21), the slope 
of the best-fit line is equal to -2 n. The intercept of the best-fit line is equal to 2/;(//;(Z, )). 

To determine the value of Z Q , equation (21) can be used again. In this case, the value of the tensile 
stress where the stress-strain curve becomes nonlinear for a particular constant strain rate tensile test is 
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used, as well as the value of the approximate inelastic strain rate when the stress-strain curve becomes 
nonlinear. The point where the stress-strain curve becomes nonlinear is defined as the approximate point 
where the curve appreciably deviates from a linear extrapolation of the initial data. The strain rate used in 
the constant strain rate test divided by 1 00 was found by trial and error to approximate the inelastic strain 
rate at this point reasonably well. Using this data, equation (21) is solved for Z, which is assumed to be 
equal to the value of Z a . Compression data can also be used, particularly for the case when the tensile 
stress at which the stress-strain curve becomes nonlinear is approximated due to the differences in the 
initial tension and compression responses as discussed above. In this case, equations (20) and (21) need to 
be reformulated based on equation (12) given that compressive stresses are being considered as opposed 
to tensile stresses. Another point to note is that using the data from the lowest strain rate test available has 
been found to give adequate values of Z Q . However, the calculations can be made using data from all the 
available strain rates, and an average taken if required to obtain the value of the constant. 

To determine the value for q for equations (13), (14), and (15), first equation (13) is integrated for the 
case of pure tensile loading (again since for woven ceramic matrix composites the tensile curves are more 
likely to display a distinct saturation stress), resulting in the following relation 

Z = Zi-(Zi-Z 0 )exp(-< ? £ e / ) (22) 

where e/ is the effective inelastic strain, computed using equation (10). At saturation, the value of the 
internal stress Z is assumed to approach Z h resulting in the exponential term approaching zero. Assuming 
that saturation occurs when the following condition is satisfied 


( 

exp 

V 


1 + V3cq 


y 


j 


0.01 


(23) 


the equation is solved for q, where e/ is the inelastic tensile strain at saturation. Also, in equation (23) the 
effective inelastic strain at saturation is simplified to the case of uniaxial tensile loading through the use 
of equations (9) and (10). If the inelastic strain at saturation is found to vary with strain rate, the 
parameter q is computed at each strain rate and regression techniques are utilized to determine an 
expression for the variation of q. If equations (14) and (15) are integrated, an expression similar to 
equation (23) is obtained. At saturation, the values of a and (3 are assumed to approach oq and (3i, and 
equations identical to equation (23) are obtained, which would lead to the same value for q. Therefore, 
identical values of q are used in equations (13), (14), and (15). 

To find the value of k, equation (10) is used along with results from a uniaxial tensile test and a pure 
shear test. Specifically, the effective inelastic strain at saturation for a tensile test at a particular total 
strain rate (usually the lowest) is set equal to the effective inelastic strain at saturation from a pure shear 
test at the same approximate total strain rate, resulting in the following expression which can be solved 
for k 


4 _ 

1 + V3cq y[?> 


(24) 


where e/ is once again the inelastic tensile strain at saturation and y/ is the inelastic engineering shear 
strain at saturation. 
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Finite Element Implementation 


The constitutive equations described above have been implemented into the commercially available 
transient dynamic finite element code LS-DYNA (ref. 10) through the use of the user defined subroutine 
(UMAT) option. For this work, the user defined routine was designed to work with shell elements, 
therefore the equations were specialized for the case of plane stress in the normal directions. A composite 
is typically a thin structure, which justifies the use of shell elements. The use of shell elements also allows 
greater flexibility in terms of defining a laminate of many composite plies, and in specifying the angle of 
orientation of the material. Furthermore, the use of shell elements provides greater flexibility and 
computational efficiency in terms of analyzing a complex structure. All of the material constants are 
entered into the input deck. The current values of the stress and any history variables at time t, along with 
the strain increments in progressing from time t to time t + At, are passed into the routine by LS-DYNA. 
The stresses and values of the history variables at time t + At are to be computed within the routine. The 
history variables for the finite element analysis are defined to include the deviatoric stresses, the 
hydrostatic stress, the values of the state variables and the inelastic strain rate tensor. To integrate the 
constitutive equations, a fourth order Runge-Kutta integrator (ref. 15) has been implemented into the 
UMAT routines. For this class of problems, implicit integration routines have often been used because of 
their inherent numerical stability (ref. 1 1). Flowever, in transient dynamic finite codes, explicit integration 
schemes are employed so an explicit integrator was required for the user defined subroutine developed in 
this work. 

In the UMAT routine, the material constants and the current values of the history variables are read 
into the routine. The values of the elastic stiffness tensor are computed based on the elastic modulus and 
Poisson’s ratio. For each step of the Runge-Kutta integration, the effective stress is computed employing 
the current stress values and the current value of the state variables a and (3 using equation (6). Next, the 
terms in the inelastic strain rate tensor are determined using equation (12). The effective inelastic strain 
rate is determined using equation (10), and then the state variable rates are determined using equations 
(13), (14), and (15). The appropriate integration formula and time step increment based on the progress 
through the Runge-Kutta integration routine is then used to compute updated values of the inelastic strain 
increments and the state variables. From the given total strain increments and the computed inelastic 
strain increments, along with the elastic constants, the values which are computed include the strain 
increments for the current step in the integration process, the stress increments and revised total stresses. 
In this process, the out of plane normal strain increments are calculated in a manner such that the plane 
stress condition is enforced, and the out of plane normal stress is set equal to zero. As is common for shell 
elements, a correction factor of 0.833 (ref. 10) is used in computing the transverse shear stresses. Revised 
values of the deviatoric and hydrostatic stresses determined from the computed stresses, and the 
integration process is continued as required. Finally, the history variables are updated, and the routine 
returns. 


Sample Calculations 

To demonstrate the ability of the developed constitutive equations to correctly simulate the nonlinear 
deformation response of woven ceramic matrix composites, a representative plain weave woven SiC/SiC 
system, with a Tyranno fiber in a SiC matrix, was analyzed. The experimental data were obtained by 
Jacobsen and Brondsted (ref. 5), in which the details of the material and the experimental procedures can 
be found. Tensile, compressive and in-plane shear stress-strain curves obtained at a quasi-static strain rate 
and room temperature are shown in figures 1, 2, and 3, respectively. The tensile and compressive results 
were measured along the fiber direction. As can be seen in the figures, the compression stress-strain curve 
is nearly linear elastic until failure, with only a slight nonlinearity. However, the tensile curve is highly 
nonlinear, due to matrix microcracking and interfacial debonding, and has a much lower failure stress 
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than the compression results. The stress-strain curve in shear also displays significant amounts of 
nonlinearity. 

To determine the material constants for the constitutive equations described in the previous section, 
the elastic and shear modulus, Poisson’s ratio, yield stress (stress where the stress-strain curve becomes 
nonlinear) in tension, compression and shear, saturation stress (stress level where the stress-strain curves 
becomes horizontal) in tension, compression and shear, and saturation strain in tension and shear were 
determined from the material data. The Poisson’s ratio data was obtained from Jacobson and Brondsted 
(ref. 5). This information is presented in table 1. Note that the saturation strain in shear is the engineering 
shear strain. Also note that since the tension and compression curves did not display a distinct saturation 
stress, curve extrapolation techniques were used to estimate the saturation stress and strain values. Using 
the information given in table 1 and the procedures described in previous sections of this paper, the 
material constants for the constitutive equations were calculated and are shown in table 2. Since strain 
rate dependent data were not available for this material, based on trends observed in other reports 
(refs. 3 and 4) a tensile saturation stress of 370 MPa was assumed at a strain rate of 0. 1 /sec. The quasi- 
static test data was assumed to be obtained at a strain rate of 5 x 1(T 5 /sec. 

Computed tension, compression and shear stress-strain curves are shown in figures 1, 2, and 3, 
respectively, along with the experimental curves for comparison. As can be seen in the figures, the 
correlation between the experimental and computed results is quite good for all three conditions. 
Qualitatively, the nonlinearity observed in the tension and shear curves is captured, and the significant 
differences in the tension and compression responses are correctly captured. Quantitatively, the tensile 
stresses are slightly under predicted in the early portion of the nonlinear range of the curve, and the shear 
stresses are slightly over predicted at higher shear strains until saturation is reached, but the differences 
are not significant. Since strain rate dependent data were not available for this material, in order to 
explore the capability of the equations to simulate a rate dependent response, a tensile curve was 
generated at a strain rate of 0. 1 /sec, and is shown in figure 4 (“high rate” in the figure) along with the 
quasi-static computed tensile curve (“low rate” in the figure) for comparison. As can be seen in the figure, 
the equations do predict a noticeable rate dependence to the tensile response, which indicates that the 
current model can capture the rate dependence in the material deformation. Future efforts will include 
obtaining strain rate dependent material data for a woven ceramic matrix composite and more thoroughly 
examining the capability of the equations to correctly simulate the observed rate dependence in the 
composite deformation. 


Conclusions 

The Bodner-Partom viscoplastic state variable model (ref. 7) has been modified in order to analyze 
the strain rate dependent, nonlinear deformation of plain weave woven ceramic matrix composites. The 
differences in the tensile and compressive deformation responses have been accounted for through the 
modification of the effective stress and effective inelastic strain definitions given in the original 
equations. The procedures used to compute the shear response have also been modified to account for the 
fact that in an orthotropic composite the shear response is independent of the normal response. The 
material constants can be obtained based on tension, compression and shear tests. The constitutive 
equations have been implemented within the transient dynamic finite element code LS-DYNA through 
the use of the user defined subroutine (UMAT) option. Sample calculations have been conducted using a 
representative material which demonstrates the ability of the model to qualitatively capture the 
nonlinearity, strain rate dependence and tension/compression asymmetry observed in the deformation 
response of woven ceramic matrix composites. 

Future efforts will be concentrated in several areas. Primarily, failure criteria will be developed and 
added to the model, in order to provide the ability to model the full deformation, damage and failure of 
woven ceramic matrix composites. In addition, as mentioned above the proper modeling of the transverse 
stresses will be considered in greater detail. 


NASA/TM— 2004-213125 


10 



References 


1. Camus, G.: “Modeling of the mechanical behavior and damage processes of fibrous ceramic matrix 
composites: application to a 2-D SiC/SiC.” International Journal of Solids and Structures, vol. 37, 
pp. 919-942, 2000. 

2. Camus, G.: Guillaumat, L.; and Baste, S.: “Development of Damage in a 2D Woven C/SiC Composite 
Under Mechanical Loading: I. Mechanical Characterization.” Composites Science and Technology, vol. 
56, pp. 1363-1372, 1996. 

3. Zhu, S.; Cao, J.-W.; Mizuno, M.; and Kagawa, Y.: “Effect of loading rate and temperature on monotonic 
tensile behavior in an enhanced SiC/SiC composite.” Scripta Materialia, vol. 50, pp. 349-352, 2004. 

4. Lipetzky, P.; Dvorak, G.J.; and Stoloff, N.S.: “Tensile properties of a SiCf/SiC composite.” Materials 
Science and Engineering, vol. A216, pp. 1 1-19, 1996. 

5. Jacobson, T.K.; and Brondsted, P.: “Mechanical Properties of Two Plain-Woven Chemical Vapor 
Infiltrated Silicon Carbide-Matrix Composites.” Journal of the American Ceramic Society, vol. 84, 
pp. 1043-1051,2001. 

6. Gasser, A.; Ladeveze, P.; and Poss, M.: “Damage Mechanisms of a Woven SiC/SiC Composite: Modeling 
and Identification.” Composites Science and Technology, vol. 56, pp. 779-784, 1996. 

7. Bodner, S.R.: Unified Plasticity for Engineering Applications . Kluwer Academic/Plenum Publishers, New 
York, 2002. 

8. Goldberg, R.K.; Roberts, G.D.; and Gilat, A: “Implementation of an Associative Flow Rule Including 
Hydrostatic Stress Effects Into the High Strain Rate Deformation Analysis of Polymer Matrix 
Composites.” NASA TM-2003-212382, National Aeronautics and Space Administration, 2003. 

9. Khan, A.S.; and Huang, S.: Continuum Theory of Plasticity . John Wiley and Sons, Inc., New York, 1995. 

10. Anonymous: LS-DYNA Keyword User’s Manual, Version 970 . Livermore Software Technology 
Corporation, Livermore, CA, 2003. 

11. Stouffer, D.C.; and Dame, L.T.: Inelastic Deformation of Metals. Models, Mechanical Properties and 
Metallurgy . John Wiley and Sons, New York, 1996. 

12. Li, F.Z; and Pan, J.: “Plane-Stress Crack-Tip Fields for Pressure-Sensitive Dilatant Materials.” Journal of 
Applied Mechanics, vol. 57, pp. 4CM9, 1990. 

13. Chang, W.J.; and Pan, J.: “Effects of Yield Surface Shape and Round-Off Vertex on Crack-Tip Fields for 
Pressure-Sensitive Materials.” International Journal of Solids and Structures, vol. 34, 

pp. 3291-3320, 1997. 

14. Hsu, S.-Y., Vogler, T.J.; and Kyriakides, S.: “Inelastic behavior of an AS4/PEEK composite under 
combined transverse compression and shear. Part II: modeling.” International Journal of Plasticity, vol. 15, 
pp. 807-836, 1999. 

15. Kreyszig, E.: Advanced Engineering Mathematics, 7 th Edition. John Wiley and Sons, Inc., New York, 
1992. 


TABLE 1.— MATERIAL DATA TO BE USED IN DETERMINING MODEL CONSTANTS 



Modulus 

(GPa) 

Poisson’s 

Ratio 

Yield Stress 
(MPa) 

Saturation 
Stress (MPa) 

Saturation 

Strain 

Tension 

168.11 

0.20 

35.00 

335.27 

0.00703 

Compression 

168.11 

0.20 

335.29 

959.38 


Shear 

64.00 


32.00 

182.00 

0.00550 


TABLE 2 . — MATERIAL CONSTANTS FOR CONSTITUTIVE EQUATIONS 


D a 

n 

Z! MPa 

Z D MPa 

7 

«i 

a 0 

Pi 

P 0 

K 

1 x 10 6 

1.91 

1374.96 

150.23 

1000.50 

0.278 

0.468 

2.484 

1.308 

3.00 
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Compressive Stress (MPa) Stress (MPa) 



Strain 

Figure 1. — Experimental and computed tensile stress-strain curves for 
woven SiC/SiC composite under quasi-static loading conditions. 



Compressive Strain 

Figure 2. — Experimental and computed compressive stress-strain curves for 
woven SiC/SiC composite under quasi-static loading conditions. 
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Stress (MPa) Shear Stress (MPa) 



Shear Strain 


Figure 3. — Experimental and computed shear stress-shear strain curves 
for woven SiC/SiC composite under quasi-static loading conditions. 



Strain 

Figure 4. — Computed tensile stress-strain curves of woven SiC/SiC composite 
at strain rates of 5 x 1(T 5 /sec (Low Rate) and 0.1 /sec (High Rate). 
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